Team 2 File
Michael Egle, John Chandara, Jaydon Cobb, Anthony Song
4/22/2020
library(tidyverse)## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.0 ──## ✓ ggplot2 3.2.1 ✓ purrr 0.3.3
## ✓ tibble 2.1.3 ✓ dplyr 0.8.4
## ✓ tidyr 1.0.2 ✓ stringr 1.4.0
## ✓ readr 1.3.1 ✓ forcats 0.5.0## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()library(lubridate)##
## Attaching package: 'lubridate'## The following object is masked from 'package:base':
##
## datelibrary(ggpubr)## Loading required package: magrittr##
## Attaching package: 'magrittr'## The following object is masked from 'package:purrr':
##
## set_names## The following object is masked from 'package:tidyr':
##
## extractlibrary(readxl)
library(ggplot2)Michael Egle Analysis
case <- read_csv(url('https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_covid19_recovered_global.csv'))## Parsed with column specification:
## cols(
## .default = col_double(),
## `Province/State` = col_character(),
## `Country/Region` = col_character()
## )## See spec(...) for full column specifications.death <- read_csv(url('https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_covid19_deaths_global.csv'))## Parsed with column specification:
## cols(
## .default = col_double(),
## `Province/State` = col_character(),
## `Country/Region` = col_character()
## )
## See spec(...) for full column specifications.rec <- read_csv(url('https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_covid19_recovered_global.csv'))## Parsed with column specification:
## cols(
## .default = col_double(),
## `Province/State` = col_character(),
## `Country/Region` = col_character()
## )
## See spec(...) for full column specifications.mob <- read_csv('applemobilitytrends-2020-04-24.csv')## Parsed with column specification:
## cols(
## .default = col_double(),
## geo_type = col_character(),
## region = col_character(),
## transportation_type = col_character()
## )
## See spec(...) for full column specifications.nyc <- read_csv('nyc_covid.csv') #from nyc.gov## Parsed with column specification:
## cols(
## DATE_OF_INTEREST = col_character(),
## Cases = col_double(),
## Hospitalizations = col_character(),
## Deaths = col_character()
## )counties <- read_csv(url('https://raw.githubusercontent.com/nytimes/covid-19-data/master/us-counties.csv'))## Parsed with column specification:
## cols(
## date = col_date(format = ""),
## county = col_character(),
## state = col_character(),
## fips = col_character(),
## cases = col_double(),
## deaths = col_double()
## )mob <- mob %>%
pivot_longer(`2020-01-13`:`2020-04-24`,
names_to = "date",
values_to = "percent_change")counties <- counties %>%
mutate(date = ymd(date))
mob <- mob %>%
mutate(date = ymd(date))Break up the county data set into four datasets for our cities of interest
nyc <- counties %>%
filter(county == "New York City") %>%
select(-fips)
mia <- counties %>%
filter(county == "Miami-Dade") %>%
select(-fips)
chi <- counties %>%
filter(county == "Cook", state == "Illinois") %>%
select(-fips)
sf <- counties %>%
filter(county == "San Francisco") %>%
select(-fips)Use inner_join and some date manipulations to get the cases from the previous day
nyc2 <- nyc
mia2 <- mia
chi2 <- chi
sf2 <- sf
nyc <- nyc %>%
mutate(day_before = date - 1)
mia <- mia %>%
mutate(day_before = date - 1)
chi <- chi %>%
mutate(day_before = date - 1)
sf <- sf %>%
mutate(day_before = date - 1)
nyc2 <- nyc2 %>%
select(c('date', 'deaths'))
mia2 <- mia2 %>%
select(c('date', 'deaths'))
chi2 <- chi2 %>%
select(c('date', 'deaths'))
sf2 <- sf2 %>%
select(c('date', 'deaths'))
colnames(nyc2) <- c("day_before", "deaths_before")
colnames(mia2) <- c("day_before", "deaths_before")
colnames(chi2) <- c("day_before", "deaths_before")
colnames(sf2) <- c("day_before", "deaths_before")
nyc <- nyc %>%
inner_join(nyc2)## Joining, by = "day_before"mia <- mia %>%
inner_join(mia2)## Joining, by = "day_before"chi <- chi %>%
inner_join(chi2)## Joining, by = "day_before"sf <- sf %>%
inner_join(sf2)## Joining, by = "day_before"nyc <- nyc %>%
mutate(new_deaths = deaths - deaths_before)
mia <- mia %>%
mutate(new_deaths = deaths - deaths_before)
chi <- chi %>%
mutate(new_deaths = deaths - deaths_before)
sf <- sf %>%
mutate(new_deaths = deaths - deaths_before)nyc_travel <- mob %>%
filter(region == "New York City")
mia_travel <- mob %>%
filter(region == "Miami")
chi_travel <- mob %>%
filter(region == "Chicago")
sf_travel <- mob %>%
filter(region == "San Francisco - Bay Area")Put together graph for NYC
nyc_travel <- nyc_travel %>%
filter(date > ymd("2020/2/29"))nyc_travel_viz <- nyc_travel %>%
ggplot(aes(x = date, y = percent_change, group = transportation_type,
color = transportation_type)) +
geom_line() +
geom_smooth() +
theme(legend.position = "none") +
annotate(geom = "text", x = ymd("2020/4/17"), y = 8, label = "Transit",
size = 3) +
annotate(geom = "text", x = ymd("2020/4/17"), y = 35, label = "Walking",
size = 3) +
annotate(geom = "text", x = ymd("2020/4/17"), y = 55, label = "Driving",
size = 3) +
geom_vline(xintercept = ymd("2020/3/20"), linetype = "dashed") +
annotate(geom = "label", x = ymd("2020/3/20"), y = 125,
label = "Stay At Home Order Put Into Effect 3/20", size = 3) +
scale_x_date(limits = c(ymd("2020-03-01", "2020-04-24"))) +
scale_y_continuous(limits = c(0, 200)) +
labs(title = "New York City",
subtitle = "26,403 people per square mile",
y = "Percent of Average Movement") +
theme(axis.title.x = element_blank())nyc_deaths_viz <- nyc %>%
ggplot(aes(x = date, y = new_deaths)) +
geom_histogram(stat = "identity") +
geom_vline(xintercept = ymd("2020/3/20"), linetype = "dashed") +
scale_x_date(limits = c(ymd("2020-03-01", "2020-04-24"))) +
geom_smooth() +
labs(caption = "Movement data from Apple Maps | Deaths data from New York Times",
x = "Date",
y = "New Deaths") +
annotate(geom = "text", x = ymd("2020/4/12"), y = 1150,
label = "11,157 total deaths as of 4/24", size = 3)## Warning: Ignoring unknown parameters: binwidth, bins, padggarrange(nyc_travel_viz, nyc_deaths_viz,
nrow = 2)## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'## Warning: Removed 8 rows containing non-finite values (stat_smooth).## Warning: Removed 8 rows containing missing values (position_stack).## Warning: Removed 1 rows containing missing values (geom_bar).
Now let’s put together Miami’s graphic
mia_travel <- mia_travel %>%
filter(date > ymd("2020/2/29"))mia_travel_viz <- mia_travel %>%
ggplot(aes(x = date, y = percent_change, group = transportation_type,
color = transportation_type)) +
geom_line() +
geom_smooth() +
theme(legend.position = "none") +
annotate(geom = "text", x = ymd("2020/4/17"), y = 25, label = "Transit",
size = 3) +
annotate(geom = "text", x = ymd("2020/4/17"), y = 40, label = "Walking",
size = 3) +
annotate(geom = "text", x = ymd("2020/4/17"), y = 55, label = "Driving",
size = 3) +
geom_vline(xintercept = ymd("2020/3/26"), linetype = "dashed") +
annotate(geom = "label", x = ymd("2020/3/26"), y = 140,
label = "Stay At Home Order Put Into Effect for Miami-Dade County 3/26",
size = 3) +
scale_x_date(limits = c(ymd("2020-03-01", "2020-04-24"))) +
scale_y_continuous(limits = c(0, 200)) +
labs(title = "Miami",
subtitle = "12,139 people per square mile",
y = "Percent of Average Movement") +
theme(axis.title.x = element_blank())mia_deaths_viz <- mia %>%
ggplot(aes(x = date, y = new_deaths)) +
geom_histogram(stat = "identity") +
geom_vline(xintercept = ymd("2020/3/26"), linetype = "dashed") +
scale_x_date(limits = c(ymd("2020-03-01", "2020-04-24"))) +
geom_smooth() +
scale_y_continuous(limits = c(0, 100)) +
labs(caption = "Movement data from Apple Maps | Deaths data from New York Times",
x = "Date",
y = "New Deaths") +
annotate(geom = "text", x = ymd("2020/4/12"), y = 70,
label = "287 total deaths as of 4/24", size = 3)## Warning: Ignoring unknown parameters: binwidth, bins, padggarrange(mia_travel_viz, mia_deaths_viz,
nrow = 2)## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'## Warning: Removed 8 rows containing non-finite values (stat_smooth).## Warning: Removed 8 rows containing missing values (position_stack).## Warning: Removed 1 rows containing missing values (geom_bar).## Warning: Removed 16 rows containing missing values (geom_smooth).
Now let’s do Chicago
chi_travel <- chi_travel %>%
filter(date > ymd("2020/2/29"))chi_travel_viz <- chi_travel %>%
ggplot(aes(x = date, y = percent_change, group = transportation_type,
color = transportation_type)) +
geom_line() +
geom_smooth() +
theme(legend.position = "none") +
annotate(geom = "text", x = ymd("2020/4/17"), y = 15, label = "Transit",
size = 3) +
annotate(geom = "text", x = ymd("2020/4/17"), y = 42, label = "Walking",
size = 3) +
annotate(geom = "text", x = ymd("2020/4/17"), y = 69, label = "Driving",
size = 3) +
geom_vline(xintercept = ymd("2020/3/20"), linetype = "dashed") +
annotate(geom = "label", x = ymd("2020/3/20"), y = 170,
label = "Stay At Home Order Put Into Effect 3/20", size = 3) +
scale_x_date(limits = c(ymd("2020-03-01", "2020-04-24"))) +
scale_y_continuous(limits = c(0, 200)) +
labs(title = "Chicago",
subtitle = "11,943 people per square mile",
y = "Percent of Average Movement") +
theme(axis.title.x = element_blank())chi_deaths_viz <- chi %>%
ggplot(aes(x = date, y = new_deaths)) +
geom_histogram(stat = "identity") +
geom_vline(xintercept = ymd("2020/3/20"), linetype = "dashed") +
scale_x_date(limits = c(ymd("2020-03-01", "2020-04-24"))) +
geom_smooth() +
labs(caption = "Movement data from Apple Maps | Deaths data from New York Times",
x = "Date",
y = "New Deaths") +
scale_y_continuous(limits = c(0, 125)) +
annotate(geom = "text", x = ymd("2020/4/12"), y = 112,
label = "1,220 total deaths as of 4/24", size = 3)## Warning: Ignoring unknown parameters: binwidth, bins, padggarrange(chi_travel_viz, chi_deaths_viz,
nrow = 2)## `geom_smooth()` using method = 'loess' and formula 'y ~ x'## Warning: Removed 1 rows containing non-finite values (stat_smooth).## `geom_smooth()` using method = 'loess' and formula 'y ~ x'## Warning: Removed 44 rows containing non-finite values (stat_smooth).## Warning: Removed 44 rows containing missing values (position_stack).## Warning: Removed 2 rows containing missing values (geom_bar).## Warning: Removed 22 rows containing missing values (geom_smooth).
And finally San Francisco
sf <- sf %>%
filter(date > ymd("2020/2/29"))sf_travel_viz <- sf_travel %>%
ggplot(aes(x = date, y = percent_change, group = transportation_type,
color = transportation_type)) +
geom_line() +
geom_smooth() +
theme(legend.position = "none") +
annotate(geom = "text", x = ymd("2020/4/17"), y = 10,
label = "Transit",
size = 3) +
annotate(geom = "text", x = ymd("2020/4/17"), y = 32,
label = "Walking",
size = 3) +
annotate(geom = "text", x = ymd("2020/4/17"), y = 52,
label = "Driving",
size = 3) +
geom_vline(xintercept = ymd("2020/3/16"), linetype = "dashed") +
annotate(geom = "label", x = ymd("2020/3/16"), y = 125,
label = "Stay At Home Order Put Into Effect 3/16",
size = 3) +
scale_x_date(limits = c(ymd("2020-03-01", "2020-04-24"))) +
scale_y_continuous(limits = c(0, 200)) +
labs(title = "San Francisco",
subtitle = "17,246 people per square mile",
y = "Percent of Average Movement") +
theme(axis.title.x = element_blank())sf_deaths_viz <- sf %>%
ggplot(aes(x = date, y = new_deaths)) +
geom_histogram(stat = "identity") +
geom_vline(xintercept = ymd("2020/3/16"), linetype = "dashed") +
scale_x_date(limits = c(ymd("2020-03-01", "2020-04-24"))) +
geom_smooth() +
scale_y_continuous(limits = c(0, 100)) +
labs(caption = "Movement data from Apple Maps | Deaths data from New York Times",
x = "Date",
y = "New Deaths") +
annotate(geom = "text", x = ymd("2020/4/12"), y = 50,
label = "22 total deaths as of 4/24", size = 3)## Warning: Ignoring unknown parameters: binwidth, bins, padggarrange(sf_travel_viz, sf_deaths_viz,
nrow = 2)## `geom_smooth()` using method = 'loess' and formula 'y ~ x'## Warning: Removed 144 rows containing non-finite values (stat_smooth).## Warning: Removed 144 rows containing missing values (geom_path).## `geom_smooth()` using method = 'loess' and formula 'y ~ x'## Warning: Removed 8 rows containing non-finite values (stat_smooth).## Warning: Removed 8 rows containing missing values (position_stack).## Warning: Removed 2 rows containing missing values (geom_bar).## Warning: Removed 14 rows containing missing values (geom_smooth).
Anthony Song Analysis
Read in the data
mob <- read_csv('applemobilitytrends-2020-04-24.csv')## Parsed with column specification:
## cols(
## .default = col_double(),
## geo_type = col_character(),
## region = col_character(),
## transportation_type = col_character()
## )## See spec(...) for full column specifications.pivot longer and month year date
mob <- mob %>%
pivot_longer(`2020-01-13`:`2020-04-24`,
names_to = "date",
values_to = "percent_change")
mob <- mob %>%
mutate(date = ymd(date))Filter some low rate countries
mob1 <- mob %>%
filter(transportation_type == "walking")
# Very low death rate per million people
indo <- mob1 %>%
filter(region == "Indonesia")
cam <- mob1 %>%
filter(region == "Cambodia")
kor <- mob1 %>%
filter(region == "Republic of Korea")
tai <- mob1 %>%
filter(region == "Taiwan")
tha <- mob1 %>%
filter(region == "Thailand")
low_rate <- rbind(cam, tha)First Visual
reg_viz3 <- low_rate%>%
ggplot(aes(x = date, y = percent_change))+geom_line()+facet_wrap(~region)+geom_smooth()+ylim(0,200)
reg_viz3+labs(x = "Date", y= "Percent Mobility Relative to Baseline", title = "Mobility of 2 Low Death Rate Regions")+geom_hline(yintercept = 125, linetype = "dashed", color = "red")+
annotate(geom = "label", x = ymd("2020/4/10"), y = 125,
label = "125% Mobility",
size = 3)## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
#These Countries have very low death rate per million, their rate small rate of change and activity agree with this.
# Rate of change of Indonesia
qplot(date,percent_change,data=indo) + stat_smooth(aes(outfit=fit_indo<<-..y..))+annotate("rect", ymin = 25, ymax = 100, xmin = ymd("2020/3/1"), xmax = ymd("2020/4/1"), alpha = .2)+labs(title = "Indonesia Rate of Change")## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_indo <- data.frame(fit_indo)
fit_indo$day<- seq.int(nrow(fit_indo))
vc_indo <- c(max(fit_indo$fit_indo), min(fit_indo$fit_indo) )
max_min_indo <- fit_indo[fit_indo$fit_indo %in% vc_indo,]
print("Indonesia ROC")## [1] "Indonesia ROC"(max_min_indo[2,1] - max_min_indo[1,1]) / (max_min_indo[2,2] -max_min_indo[1,2])## [1] -1.451253# rate of change of Cambodia
qplot(date,percent_change,data=cam) + stat_smooth(aes(outfit=fit_cam<<-..y..))## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_cam <- data.frame(fit_cam)
fit_cam$day<- seq.int(nrow(fit_cam))
vc_cam <- c(max(fit_cam$fit_cam), min(fit_cam$fit_cam) )
max_min_cam <- fit_cam[fit_cam$fit_cam %in% vc_cam,]
print("Camboia ROC")## [1] "Camboia ROC"(max_min_cam[2,1] - max_min_cam[1,1]) / (max_min_cam[2,2] -max_min_cam[1,2])## [1] -0.7241492# rate of change of Soutn korea
qplot(date,percent_change,data=kor) + stat_smooth(aes(outfit=fit_kor<<-..y..))## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_kor <- data.frame(fit_kor)
fit_kor$day<- seq.int(nrow(fit_kor))
vc_kor <- c(max(fit_kor$fit_kor), min(fit_kor$fit_kor) )
max_min_kor <- fit_kor[fit_kor$fit_kor %in% vc_kor,]
print("South Korea ROC")## [1] "South Korea ROC"(max_min_kor[2,1] - max_min_kor[1,1]) / (max_min_kor[2,2] -max_min_kor[1,2])## [1] -1.500561# rate of change of Taiwan
qplot(date,percent_change,data=tai) + stat_smooth(aes(outfit=fit_tai<<-..y..))## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_tai <- data.frame(fit_tai)
fit_tai$day<- seq.int(nrow(fit_tai))
vc_tai <- c(max(fit_tai$fit_tai), min(fit_tai$fit_tai) )
max_min_tai <- fit_tai[fit_tai$fit_tai %in% vc_tai,]
print("Taiwan ROC")## [1] "Taiwan ROC"(max_min_tai[2,1] - max_min_tai[1,1]) / (max_min_tai[2,2] -max_min_tai[1,2])## [1] -0.5122063# rate of change of Thailand
qplot(date,percent_change,data=tha) + stat_smooth(aes(outfit=fit_tha<<-..y..))## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_tha <- data.frame(fit_tha)
fit_tha$day<- seq.int(nrow(fit_tha))
vc_tha <- c(max(fit_tha$fit_tha), min(fit_tha$fit_tha) )
max_min_tha <- fit_tha[fit_tha$fit_tha %in% vc_tha,]
print("Thailand ROC")## [1] "Thailand ROC"(max_min_tha[2,1] - max_min_tha[1,1]) / (max_min_tha[2,2] -max_min_tha[1,2])## [1] -1.275496Regions
#High death rate per million people.
ity <- mob %>%
filter(region == "Italy")
por <- mob %>%
filter(region == "UK")
esp <- mob %>%
filter(region == "Spain")
bel <- mob %>%
filter(region == "Belgium")
ire <- mob %>%
filter(region == "Ireland")
usa <- mob %>%
filter(region == "United States")
non_rec_regions <- rbind( por, esp, ire, bel)
non_rec_regions <- non_rec_regions %>%
filter(transportation_type == "walking")
reg_viz2 <- non_rec_regions %>%
ggplot(aes(x = date, y = percent_change))+geom_line()+facet_wrap(~region)+geom_smooth()+ylim(0,200)
reg_viz2+labs(x = "Date", y = "Percent Change in Mobility", title = "Mobility of 4 High Death Rate Regions")+labs(x = "Date", y= "Percent Change in Mobility")+geom_hline(yintercept = 125, linetype = "dashed", color = "red")+annotate(geom = "label", x = ymd("2020/4/10"), y = 125,
label = "125% Mobility",
size = 3)## `geom_smooth()` using method = 'loess' and formula 'y ~ x'## Warning: Removed 3 rows containing non-finite values (stat_smooth).
# How does the USA compare???
reg_viz_USA <- usa %>%
ggplot(aes(x = date, y = percent_change))+geom_line()+facet_wrap(~region)+geom_smooth()+ylim(0,175)
reg_viz_USA+labs(x = "Date", y = "Percent Change in Mobility", title = "Walking Mobility of the USA")+labs(x = "Date", y= "Percent Mobility Relative to Baseline")+geom_hline(yintercept = 125, linetype = "dashed", color = "red")+annotate(geom = "label", x = ymd("2020/4/10"), y = 125,
label = "125% Mobility",
size = 3)## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
high_low <-rbind(cam, tha, bel, ire)
reg_viz5 <- high_low%>%
ggplot(aes(x = date, y = percent_change))+geom_line()+facet_wrap(~region)+geom_smooth()+ylim(0,200)+labs(x = "Date", y= "Percent Mobility Relative to Baseline", title = "Contrasting Walking Mobility Data of 4 Regions")+geom_hline(yintercept = 125, linetype = "dashed", color = "red")+ annotate(geom = "label", x = ymd("2020/4/10"), y = 125,label = "125% Mobility",size = 3)
reg_viz5## `geom_smooth()` using method = 'loess' and formula 'y ~ x'## Warning: Removed 3 rows containing non-finite values (stat_smooth).
# rate of change of italy
qplot(date,percent_change,data=ity) + stat_smooth(aes(outfit=fit_ity<<-..y..))## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_ity <- data.frame(fit_ity)
fit_ity$day<- seq.int(nrow(fit_ity))
vc_ity <- c(max(fit_ity$fit_ity), min(fit_ity$fit_ity) )
max_min_ity <- fit_ity[fit_ity$fit_ity %in% vc_ity,]
print("Italy ROC")## [1] "Italy ROC"rate <- c((max_min_ity[2,1] - max_min_ity[1,1]) / (max_min_ity[2,2] -max_min_ity[1,2]))
# rate of change portugal
qplot(date,percent_change,data=por) + stat_smooth(aes(outfit=fit_por<<-..y..))## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_por <- data.frame(fit_por)
fit_por$day<- seq.int(nrow(fit_por))
vc_por <- c(max(fit_por$fit_por), min(fit_por$fit_por) )
max_min_por <- fit_por[fit_por$fit_por %in% vc_por,]
print("Portugal ROC")## [1] "Portugal ROC"rate <- append(rate, (max_min_por[2,1] - max_min_por[1,1]) / (max_min_por[2,2] -max_min_por[1,2]))
# rate of change spain
qplot(date,percent_change,data=esp) + stat_smooth(aes(outfit=fit_esp<<-..y..))## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_esp <- data.frame(fit_esp)
fit_esp$day<- seq.int(nrow(fit_esp))
vc_esp <- c(max(fit_esp$fit_esp), min(fit_esp$fit_esp) )
max_min_esp <- fit_esp[fit_esp$fit_esp %in% vc_esp,]
print("Spain ROC")## [1] "Spain ROC"rate <- append(rate,(max_min_esp[2,1] - max_min_esp[1,1]) / (max_min_esp[2,2] -max_min_esp[1,2]))
# rate of change ireland
qplot(date,percent_change,data=ire) + stat_smooth(aes(outfit=fit_ire<<-..y..))## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_ire <- data.frame(fit_ire)
fit_ire$day<- seq.int(nrow(fit_ire))
vc_ire <- c(max(fit_ire$fit_ire), min(fit_ire$fit_ire) )
max_min_ire <- fit_ire[fit_ire$fit_ire %in% vc_ire,]
print("Ireland ROC")## [1] "Ireland ROC"rate <- append(rate,(max_min_ire[2,1] - max_min_ire[1,1]) / (max_min_ire[2,2] -max_min_ire[1,2]))
# rate of change Belgium
qplot(date,percent_change,data=bel) + stat_smooth(aes(outfit=fit_bel<<-..y..))## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_bel <- data.frame(fit_bel)
fit_bel$day<- seq.int(nrow(fit_bel))
vc_bel <- c(max(fit_bel$fit_bel), min(fit_bel$fit_bel) )
max_min_bel <- fit_bel[fit_bel$fit_bel %in% vc_bel,]
print("Belgium ROC")## [1] "Belgium ROC"rate <- append(rate,(max_min_bel[2,1] - max_min_bel[1,1]) / (max_min_bel[2,2] -max_min_bel[1,2]))
rate## [1] -2.600079 -2.151073 -3.501693 -2.762939 -2.659418# rate of change of USA
qplot(date,percent_change,data=usa) + stat_smooth(aes(outfit=fit_usa<<-..y..))## Warning: Ignoring unknown aesthetics: outfit## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
fit_usa <- data.frame(fit_usa)
fit_usa$day<- seq.int(nrow(fit_usa))
vc_usa <- c(max(fit_usa$fit_usa), min(fit_usa$fit_usa) )
max_min_usa <- fit_usa[fit_usa$fit_usa %in% vc_usa,]
max_min_usa## fit_usa day
## 29 114.7899 29
## 70 40.9266 70(max_min_usa[2,1] - max_min_usa[1,1]) / (max_min_usa[2,2] -max_min_usa[1,2])## [1] -1.801544Read excel
dpp <- read_excel("deaths_per_million.xlsx")## New names:
## * `` -> ...2
## * `` -> ...3
## * `` -> ...4clean data
colnames(dpp) <- dpp[1, ]
dpp <- dpp[-1, ]
dpp[] <- lapply(dpp, function(x) type.convert(as.character(x)))
dpp <- dpp[complete.cases(dpp),]
# Take absolute value for plotting purposes.
#dpp$`Rate of decrease in Mobility` <- abs(dpp$`Rate of decrease in Mobility`)
head(dpp)## # A tibble: 6 x 4
## Country `Deaths per a million peop… `Rate of decrease in Mobil… `below 125%`
## <fct> <dbl> <dbl> <int>
## 1 Indonesia 3 -1.45 1
## 2 Cambodia 0 -0.724 1
## 3 Korea 5 -1.50 1
## 4 Taiwan 0.3 -0.512 1
## 5 Thailand 0.71 -1.28 1
## 6 Italy 436 -2.60 0Run linear regression
# Model a polynomial function.
regressor = lm(dpp$`Deaths per a million people` ~ poly(dpp$`Rate of decrease in Mobility`,2, raw = TRUE), data=dpp)
summary(regressor)##
## Call:
## lm(formula = dpp$`Deaths per a million people` ~ poly(dpp$`Rate of decrease in Mobility`,
## 2, raw = TRUE), data = dpp)
##
## Residuals:
## Min 1Q Median 3Q Max
## -140.399 -90.931 5.949 63.518 215.486
##
## Coefficients:
## Estimate Std. Error
## (Intercept) -127.82 161.52
## poly(dpp$`Rate of decrease in Mobility`, 2, raw = TRUE)1 -113.78 181.38
## poly(dpp$`Rate of decrease in Mobility`, 2, raw = TRUE)2 25.54 45.36
## t value Pr(>|t|)
## (Intercept) -0.791 0.452
## poly(dpp$`Rate of decrease in Mobility`, 2, raw = TRUE)1 -0.627 0.548
## poly(dpp$`Rate of decrease in Mobility`, 2, raw = TRUE)2 0.563 0.589
##
## Residual standard error: 119.1 on 8 degrees of freedom
## Multiple R-squared: 0.7754, Adjusted R-squared: 0.7193
## F-statistic: 13.81 on 2 and 8 DF, p-value: 0.002543pred_dpp <- data.frame(death = predict(regressor, dpp), rate=dpp$`Rate of decrease in Mobility`)
#Graph with R^2 and p-value
ggplot(data = dpp, aes(y = `Deaths per a million people` , x = `Rate of decrease in Mobility`)) +
geom_point(color='blue') +
geom_line(color='red',data = pred_dpp, aes(x=rate, y=death))+ labs(x = "Rate of Mobility from Peak to Base (Percent/Day)", title = "Effect of Walking Mobility Rate on COVID-19 Death Rate")+ annotate(geom = "label", x = -1, y = 500,label = "R^2: 0.775, p-value: 0.003",size = 5)+xlim(-4,0)+annotate(geom = "label", x = -1.6, y = 170,label = "USA",size = 4)+annotate(geom = "label", x = -0.8, y = 50,label = "Taiwan",size = 4)+annotate(geom = "label", x = -3.5, y = 450,label = "Spain",size = 4)
AS2*